Site description

We performed this study at the University of Michigan Biological Station (UMBS; Ameriflux site US-UMd, 45.5625\(^{\circ}\), -84.6975\(^{\circ}\)), located in Northern Lower Michigan, USA. The area surrounding the research station is 87% well-drained upland forest and 13% wetland (Bergen and Dronova 2007); the focus of this study is on the former. The landscape geography of the UMBS upland forest is 20.4% moraine, 37.8% high outwash plain, 31.3% low outwash plain, 5.7% lake-margin terrace, 3.6% ice-contact, and 1.2% lowland glacial lake (Bergen and Dronova 2007). Most of the UMBS upland forest canopy is dominated by temperate deciduous early-successional species, most importantly Populus grandidentata (bigtooth aspen) and, to a lesser extent, Betula papyrifera (paper birch), with Acer saccharum (sugar maple), Acer rubrum (red maple), Fagus grandifolia (American beech), Tilia americana (basswood), Betula alleghaniensis (yellow birch), Fraxinus americana (white ash), Tsuga canadensis (eastern hemlock), Quercus rubra (northern red oak), Pinus strobus (white pine), and Pinus resinosa (red pine) existing in various fractions in the understory (and, in patches, in the canopy). This composition is a legacy of the site’s disturbance history: The site was intensively logged in the late 1800s and early 1900s, and experienced regularly recurring fires until the mid-1920s, at which point a regime of active fire suppression started that has persisted to the present day. As a result, the average stand age in 2013 was 95 years. The majority of the forest that is aspen-dominated is undergoing succession to “northern hardwood” (maple, beech, basswood, birch, ash, hemlock), “upland conifer” (red and white pine), or “northern red-oak” ecotypes (Bergen and Dronova 2007).

Forest stands at UMBS were previously exposed to experimental disturbance as part of the Forest Accelerated Succession Experiment (FASET).

Ecosystem Demography Model (version 2.2)

For this study, we used the Ecosystem Demography Model, version 2.2 (ED2). Below, we provide only a high-level description of key model features and focus in detail only on submodels and configurations specific to this study. A full description of the default configuration of the model is provided by Longo et al. (2019b).

The key feature of the “ecosystem demography” approach underlying the ED2 model is a system of partial differential equations that simulate the emergent mean behavior of “cohorts” of individuals of similar composition, size, and age (Moorcroft et al. 2001). The original “Ecosystem Demography” (ED) model (Moorcroft et al. 2001) was enhanced by Medvigy et al. (2009) (as ED2) with fully-coupled leaf physiology, canopy radiative transfer, and micrometeorology, and each of these components have been further refined over the last decade (culminating in ED2.2; Longo et al. 2019b). Various versions of the ED model have been validated in boreal (Medvigy and Moorcroft 2011), temperate (Medvigy et al. 2009, Dietze et al. 2014, Raczka et al. 2018), and tropical biomes (Moorcroft et al. 2001, Longo et al. 2019a), and have been applied in:

• paleoecological studies (Rollinson et al. 2017);
• free-air CO₂ enrichment studies (De Kauwe et al. 2013, Miller et al. 2015);
• estimation of potential carbon stocks (Hurtt et al. 2004)
• vegetation-climate feedbacks {Swann, etc}

ED2 represents both spatial and temporal variability hierarchically. The largest spatial unit is the “polygon”, within which all elements share the same meteorological forcings. Polygons are divided into multiple “sites”, each site with its own set of abiotic characteristics such as soil texture and topography. Sites are divided into discrete “patches” according to time since last disturbance. Finally, each patch is composed of an arbitrary (and dynamic) number of “cohorts”, representing the mean behavior of groups of plant individuals of similar size and composition. Temporal variability in ED2 is similarly hierarchical. Thermodynamics, including energy, water, and eddy (e.g. CO₂) fluxes, are resolved on the order of seconds. Biophysics, including photosynthesis, respiration, and radiative transfer, and the evaluation of energy, water, and CO2 budgets, are resolved on the order of minutes. Phenology and plant carbon balance are evaluated daily. Cohort dynamics, including structural growth, mortality, and reproduction, are evaluated monthly. Finally, patch dynamics (e.g. annual disturbances) are evaluated annually. In this study, we consider only a single patch (and, therefore, a single site and polygon) starting from a “near bare ground” condition.

{Brief summary of thermodynamics?} {Brief summary of water cycle?}

ED2 represents the carbon flux for each cohort as the net sum of contributions from turbulent mixing, leaf photosynthesis and respiration, (fine) root respiration, growth and maintenance respiration, and heterotrophic respiration. In addition, each cohort also has a virtual carbon balance that links fast processes like photosynthesis and respiration to slow processes like allocation (to growth, storage, and reproduction) and mortality.

The official ED2 source code is available at https://github.com/EDmodel/ED2, and the exact version used in this study (which includes minor revisions to accommodate the requirements of this study) can be obtained https://github.com/ashiklom/ED2/tree/b048950971e91699b78fc566df133caf977fda32.

ED2 submodels

In this study, we ran a factorial combination of the following ED2 configurations: (1) two-stream vs. multiple-scatter canopy radiative transfer models; (2) infinite vs. finite crown area (a.k.a. complete vs. partial shading); and (3) static vs. light-plastic traits.

Both canopy radiative transfer models in ED2 resolve the full vertical radiation profile within a patch as a function of canopy structure (leaf and wood area indices, crown area, leaf angle distribution) and incident solar radiation. The two models differ in the mathematical approximations of the underlying differential radiative transfer equations.

• Two-stream model (ED default)
  • Same general approach as in CLM 4.5 (Oleson et al. 2013), but rather than solving for just two layers (sunlit and shaded), solve for each canopy layer.
  • Set up the two-stream equations for each layer as a linear system and solve the entire set of equations numerically
  • Similar to approaches used in other models: JULES (Mercado et al. 2007), ORCHIDEE-CAN (Naudts et al. 2015), CLM(SPA) (Bonan et al. 2014)
  • For details, see Supplements S9 and S10 in (Longo et al. 2019b).
• Multiple scatter model
  • Based on Zhao and Qualls (2005)

The crown area submodel in ED2 determines the nature of competition for light between cohorts. In the default configuration (“infinite crown area”, or “complete shading”), the leaf area of a cohort is distributed across the entire horizontal area of a patch. This means that taller cohorts always receive more incoming radiation than shorter cohorts, even when the height difference is small. Put another way, this means that there is no horizontal competition, only vertical. This has been shown to excessively suppress competition from sub-dominant individuals and result in unrealistically homogeneous canopies (Fisher et al. 2015). In the alternate configuration (“finite crown area”, or “partial shading”), canopies take up only a fraction of the available horizontal area, meaning that multiple cohorts of similar height can receive the same level of light. The horizontal area of crowns is determined by allometric equations from Dietze and Clark (2008).

The third submodel we evaluated was trait plasticity. In the default configuration, all cohorts of a given plant functional type will have the same parameters, regardless of environmental conditions. This ignores the globally-documented intraspecific trait variability as a function of light level (Niinemets 2010, Keenan and Niinemets 2016). In the alternate configuration, as light level decreases (trees become more shaded), specific leaf area increases and V_cmax decreases, following empirical relationships from the tropics (Lloyd et al. 2010).

These sub-models are summarized in Table X.

ED2 parameters

The full list of parameters in ED2 is large, with over 100 parameters per plant functional type. A full sensitivity and uncertainty analysis across all of these parameters is outside the scope of this study. Instead, we selected a subset of parameters that were identified as important by previous ED2 sensitivity studies (Dietze et al. 2014, Raczka et al. 2018, Viskari et al. 2019).

Three parameters are related to leaf-level physiology: Following the enzyme-kinetic model of Farquhar et al. (1980), the rate of photosynthesis is the minimum of light-limited and enzyme-limited reactions. The former are controlled by the quantum efficiency parameter—maximum efficiency with which absorbed photosynthetically active radiation is converted to CO₂. The latter are controlled by V_cmax, the maximum rate of carbon fixation by Rubisco. The water demand of photosynthesis is controlled by the stomatal slope, the sensitivity of stomatal conductance of CO₂ as a function of CO₂ concentration and humidity at the leaf surface (Ball et al. 1987, Leuning 1995). Three more parameters correspond to the respiration rates of leaves, roots, and “growth maintenance”. Two more parameters control carbon allocation: One is the ratio of fine root to leaf biomass (fineroot2leaf, or q), and another is the ratio of “storage” carbon allocated to reproduction (r_fract). Three parameters control various aspects of adult tree mortality, namely the density-independent mortality rate (mort3), and the time scale (mort1) and critical carbon balance (mort2) of mortality from carbon starvation. An additional parameter controls seedling mortality. Specific leaf area (SLA) is used to convert leaf biomass to leaf area index, which in turn is used in a variety of calculations related to canopy radiative transfer and micrometeorology. Two additional parameters describe canopy structure in ways critical to the canopy radiative transfer: Canopy clumping factor describes how evenly leaf area is distributed in horizontal space (1 being perfectly evenly; 0 being a “black hole” where all leaves are concentrated in a single point); and leaf orientation factor describes the average distribution of leaf angles (-1 being perfectly vertical, 1 being perfectly horizontal, and 0 being random). Four parameters control leaf optical properties, namely the fractions of light reflected or transmitted in visible and near-infrared wavelengths (leaf_(reflect|trans)_(vis|nir)). Three other parameters:

• Fraction of litter carbon that enters the labile (fast) soil carbon pool (f_labile)
• Minimum plant height for reproduction (minimum_height)
• Soil-plant water conductance, or “water availability factor” (water_conductance)

By default, ED2 supports 17 different plant functional types, which divide plant species according to photosynthetic pathway (C3 vs. C4), growth form (grass vs. tree), leaf phenology habit (deciduous vs. evergreen), biome (e.g. temperate vs. tropical), and successional status (e.g. early, mid, late). However, we limited our simulations to the four plant functional types that have any appreciable presence at UMBS: Early, mid, and late temperate deciduous trees and pines. The species comprising these plant functional types are shown in Table X.

For each plant functional type, we generated a distribution of parameter values via a trait-meta analysis (LeBauer et al. 2013, see also Dietze et al. 2014, Raczka et al. 2018). Prior distributions for this meta-analysis were the same as those used in Dietze et al. (2014), Raczka et al. (2018), and Viskari et al. (2019). Species trait data for this meta-analysis came from existing records in the BETY database (www.betydb.org, LeBauer et al. 2017), as well as from publicly available records in the TRY database (www.try-db.org, Kattge et al. 2011) and, specifically for leaf optical properties, from Shiklomanov (2019 {dissertation, chapter 3}). The resulting parameter distributions are shown in Figure X. All parameters not described in this section were set to their ED2 PFT-specific defaults.

Figure 1: Input parameter distributions from PEcAn trait meta-analysis.

Modeling workflow

For each factorial combination of model configurations/submodels (described above), we ran 100 ensemble members from 1901 to 2000. Each ensemble member was initialized from ED2’s “near-bare ground” condition: An equal number of seedlings of each plant functional type (see previous section) at the minimum resolvable size. Driving meteorological data was 6-hourly CRU-NCEP combined with an annual atmospheric CO₂ record from Law-Dome ice core (Etheridge et al. 1998) and Mauna Loa observatory (Thoning et al. 1989). Soil texture was set to 92% sand, 7% silt, and 1% clay, per site-level observations in (Gough et al. 2010). The initial soil moisture profile was set to the average soil moisture profile reported in the UMBS Ameriflux ancillary data (https://ameriflux.lbl.gov/sites/siteinfo/US-UMd). The entire model execution workflow was conducted using the Predictive Ecosystem Analyzer (PEcAn, www.pecanproject.org).

Analysis of results

• Overall summary:
  • Overall parameter uncertainty – Time-series plots of total annual GPP and NPP, and maximum LAI and AGB, averaged across ensemble members
  • Impact of parameter uncertainty on ecological trajectories – plot of max annual LAI by PFT
  • Both compared against observations of recently observed LAI and NPP at UMBS (e.g. FASET, FoRTE)
• Uncertainty analysis:
  • Parameter vs. structural uncertainty – Variance in mean annual GPP from 19XX to 19YY across ensembles within model structure, and compare against
  • Breakdown of parameter uncertainty – PEcAn sensitivity and uncertainty analysis (LeBauer et al. 2013, Dietze et al. 2014)
• Details of structural uncertainty:
  • Vertical profiles of radiation, traits?

Summary

Figure 2: ED2 plot-level ensemble predictions of gross (GPP) and net (NPP) primary productivity, total leaf area index (LAI), and aboveground biomass (AGB) by model configuration. Individual ensemble members are shown as gray lines. The solid line is a LOESS regression (ggplot2::stat_smooth). Shading represents the middle 80% (10th to 90th percentile) of predictions. Black dot with error bars is the mean and min/max value from Hardiman et al. (2013).

Model projections of gross and net primary productivity, aboveground biomass, and leaf area index varied substantially with both parameter value and model structure (Figure 2). For every model structural configuration, there were multiple ensemble members that experienced total ecosystem collapse before the end of the simulation. However, model runs with finite canopy radius and two-stream radiative transfer had the highest failure rate, with most ensemble members collapsing almost immediately and the remaining few members showing relatively low productivity that declined to nearly zero. The remaining model configurations all followed a similar trend of rapid growth in the first simulation decade followed by gradual decline through the rest of the simulation (with individual ensemble members sometimes undergoing sudden collapse in specific years). The ensemble means of these configurations had very similar ecosystem states at peak productivity and at the end of the simulation, as well as the pace of initial rapid growth. That being said, ensemble means from model configurations that included finite canopy radius and multiple-scatter radiative transfer sustained peak aboveground biomass and leaf area indices for longer than other configurations, though this was not true of gross or mean primary productivity.

For most model configurations, the ensemble distribution of net primary productivity at the end of the run overlapped observations by Hardiman et al. (2013), with an ensemble mean that was slightly lower. (The exceptions were model configurations that included the combination of finite canopy radius and two-stream radiative transfer, which consistently underestimated productivity). On the other hand, ensemble predictions of leaf area index were consistently lower than observations, with only a subset of individual ensemble members across all configurations even overlapping the lower bounds of the observations. In fact, most ensemble members across all configurations failed to reach the observed leaf area index even at their peak. A subset of specific parameter combinations from model configurations that included finite canopy radius and multiple-scatter radiative transfer accurately reproduced the mean observed leaf area index.

Figure 3: ED2 PFT-level predictions of leaf area index (LAI) by model configuration and parameterization.

Different parameter combinations led to model outputs that had not only different bulk ecosystem properties but also qualitatively different plant communities (Figure 3). Most model configurations could produce either a forest that was persistently and completely dominated by early hardwood or pine, or a forest with some coexistence of multiple PFTs. The exceptions, again, were model configurations combining finite canopy radius and the two-stream radiative transfer; in these configurations, the only outcomes were immediate ecosystem collapse of all PFTs or survival and mono-dominance of pine.

Uncertainty analysis

Comparison of parameteric uncertainty within ensembles (colored bars, colored by model type) and “structural” uncertainty (variance in ensemble means; black bar) by output variable. Output variables are expressed as growing season averages for all years after 1910.

Take-away: Overall, structural uncertainty (variance in ensemble means across structures) is comparable or lower than parameter uncertainty (uncertainty within ensembles for a given structure). Parameter uncertainty varies significantly across model structures. Enabling trait plasticity increases variance in predictions of AGB. Enabling the finite canopy radius increases the variance in LAI predictions. GPP and NPP have the highest variance by far for the combination of finite canopy radius and two-stream radiative transfer.

PEcAn-liked parameter sensitivity and uncertainty analysis, by model type. Elasticity (a) is the normalized sensitivity of the model to a fixed change in the parameter. pvar (b) is the partial variance, which describes the overall contribution of the parameter to model predictive uncertainty based on the combination of parameter uncertainty and model sensitivity. Input parameter distributions are shown in Figure XXX.

Take-away: Model parameter uncertainty is sensitive to model configuration. For example, turning on the finite canopy radius with the two-stream canopy RTM increases sensitivity across many parameters.

Details about model structures

(TODO)

pagebreak

Radiative transfer: Definitions

• \(\omega\) – Clumping factor
• \(LAI\) – Leaf area index
• TAI – Total area index (leaf and woody area)
• \(LAI_e\) / \(TAI_e\) – “Effective” area indices (accounting for clumping)

Two-stream RTM

Multiple scatter RTM

\[ LAI_e = \omega LAI \]

\[ TAI_e = LAI_e + WAI \]

\[ TAI_{l} = \frac{TAI_e}{CAI} \]

Leaf (\(w_l\)) and wood (\(w_w\)) weights:

\[ w_l = \frac{LAI_e}{TAI_e} \] \[ w_w = 1 - w_l \]

Projected area, based on coefficients \(\phi_1\) and \(\phi_2\).

\[ a_{proj} = \phi_1 + \phi_2 \mu \]

Introduction

• Why we need models?
  • Forests are important, for C sequestration, other ecosystem services
  • Policy decisions related to forests depend on future projections
  • We need models to make future projections
• What kinds of models are out there?
  • Traditional: “Big leaf” models
    • Needed to be simple to run on old computers
    • Don’t need as much data
    • …but they fail to capture important processes, which matter not only for accurate predictions about C cycle, but also because they are processes we care about (habitat, forestry, tourism, etc.)!
  • Individual-based models
    • Much more detailed, but harder to run, and need much more data.
    • Trade-off between complexity and accuracy.
• Cohort-based models like ED2 as a middle ground?
  • Core ED concept (Moorcroft et al. 2001)
  • ED2 developments – (Medvigy et al. 2009, Medvigy and Moorcroft 2011)
  • Most recent version – (Longo et al. 2019b)
  • Validation/applications of the ED model:
    • (Rollinson et al. 2017) – ED2 paleoclimate
    • (Miller et al. 2015) – ED2 at Duke FACE
    • (De Kauwe et al. 2013) – ED2 at Duke and Oak Ridge FACE
    • (Longo et al. 2019a) – Model evaluation in Tropics
• How much complexity do we need to accurately model ecological dynamics?
  • Lots of existing work on parameteric uncertainty
    • (Dietze et al. 2014)
    • (Raczka et al. 2018)
    • (Pandit et al. 2018)
  • Less on structural uncertainty – i.e. not just what coefficients do we use, but what equations, and more generally, what processes do we even consider?
• What are the processes that drive succession of ecological communities?
  • Competition for light (Bazzaz 1979)
  • Others? Leaning towards omitting because factorial combination becomes too large.
    • Competition for water?
    • Competition for nutrients? (ED’s N model isn’t very good. Alternatives are being developed, but not in ED mainline yet.)
• Brief review of community/ecosystem ecology of the Upper Midwest
  • UMBS is a useful site for testing these theories because of its natural history.
  • Around 1900, mature forest dominated by white (Pinus strobus) and red pine (Pinus resinosa), hemlock (Tsuga canadensis), and sugar maple (Acer saccharum).
  • But, repeat harvest and burning between 1880 and 1920 (Kilburn 1960) led to replacement by early successional aspen (Populus grandidentata and Populus tremuloides), birch (Betula papyrifera), and cherry (Prunus pensylvatica) (???, Gates 1930, Friedman and Reich 2005, Bergen and Dronova 2007).
  • Such aspen-dominated secondary forests are sustained by background levels of disturbance, particularly fire (Gates 1930)
  • However, due to aggressive fire suppression in the US, these forests are now undergoing succession, primarily to mid/late-successional hardwood (typified by sugar maple and American beech, Fagus grandifolia) or upland conifer (typified by red and white pine) (Bergen and Dronova 2007).
  • The combination of stand-replacing disturbance ~100 years ago and suppressed disturbance since then makes this an ideal site for testing a vegetation model run from bare ground without having to prescribe disturbance.
• Guiding questions:
  • Which of these processes do we need to be able to accurately model what happens (“community succession and C cycle dynamics”) to Midwestern forests?
  • What is the cost of considering these processes, in terms of additional parameteric uncertainty?
  • What are the relative contributions of parameteric uncertainty (data limitation) vs. structural uncertainty (theoretical limitation?)?
• Study overview:
  • ED with factorial combination of submodels related to radiative transfer formulation (two-stream vs. multiple scatter), horizontal competition (finite canopy radius vs. complete shading), and trait plasticity (whether or not SLA and Vcmax vary with light level)
  • Ensemble analysis to look at parameter uncertainty

\[ \lambda = \frac{a_{proj}}{\mu} \]

Methods

#> ─ Session info ──────────────────────────────────────────────────────────
#>  setting  value                       
#>  version  R version 3.6.0 (2019-04-26)
#>  os       macOS High Sierra 10.13.6   
#>  system   x86_64, darwin17.7.0        
#>  ui       unknown                     
#>  language (EN)                        
#>  collate  en_US.UTF-8                 
#>  ctype    en_US.UTF-8                 
#>  tz       America/New_York            
#>  date     2019-07-01                  
#> 
#> ─ Packages ──────────────────────────────────────────────────────────────
#>  package       * version    date       lib source        
#>  assertthat      0.2.1      2019-03-21 [1] CRAN (R 3.6.0)
#>  backports       1.1.4      2019-04-10 [1] CRAN (R 3.6.0)
#>  base64url       1.4        2018-05-14 [1] CRAN (R 3.6.0)
#>  callr           3.2.0      2019-03-15 [1] CRAN (R 3.6.0)
#>  cli             1.1.0      2019-03-19 [1] CRAN (R 3.6.0)
#>  codetools       0.2-16     2018-12-24 [2] CRAN (R 3.6.0)
#>  colorspace      1.4-1      2019-03-18 [1] CRAN (R 3.6.0)
#>  cowplot       * 0.9.4      2019-01-08 [1] CRAN (R 3.6.0)
#>  crayon          1.3.4      2017-09-16 [1] CRAN (R 3.6.0)
#>  desc            1.2.0      2018-05-01 [1] CRAN (R 3.6.0)
#>  devtools        2.0.2      2019-04-08 [1] CRAN (R 3.6.0)
#>  digest          0.6.19     2019-05-20 [1] CRAN (R 3.6.0)
#>  dplyr         * 0.8.1      2019-05-14 [1] CRAN (R 3.6.0)
#>  drake         * 7.4.0      2019-06-07 [1] CRAN (R 3.6.0)
#>  evaluate        0.14       2019-05-28 [1] CRAN (R 3.6.0)
#>  forcats       * 0.4.0      2019-02-17 [1] CRAN (R 3.6.0)
#>  fortebaseline * 0.0.0.9000 2019-06-24 [1] local         
#>  fs            * 1.3.1      2019-05-06 [1] CRAN (R 3.6.0)
#>  fst           * 0.9.0      2019-04-09 [1] CRAN (R 3.6.0)
#>  furrr         * 0.1.0      2018-05-16 [1] CRAN (R 3.6.0)
#>  future        * 1.13.0     2019-05-08 [1] CRAN (R 3.6.0)
#>  future.callr  * 0.4.0      2019-01-07 [1] CRAN (R 3.6.0)
#>  ggplot2       * 3.2.0      2019-06-16 [1] CRAN (R 3.6.0)
#>  globals         0.12.4     2018-10-11 [1] CRAN (R 3.6.0)
#>  glue            1.3.1      2019-03-12 [1] CRAN (R 3.6.0)
#>  gtable          0.3.0      2019-03-25 [1] CRAN (R 3.6.0)
#>  here          * 0.1        2017-05-28 [1] CRAN (R 3.6.0)
#>  highr           0.8        2019-03-20 [1] CRAN (R 3.6.0)
#>  hms             0.4.2      2018-03-10 [1] CRAN (R 3.6.0)
#>  htmltools       0.3.6      2017-04-28 [1] CRAN (R 3.6.0)
#>  igraph          1.2.4.1    2019-04-22 [1] CRAN (R 3.6.0)
#>  knitr           1.23       2019-05-18 [1] CRAN (R 3.6.0)
#>  labeling        0.3        2014-08-23 [1] CRAN (R 3.6.0)
#>  lattice         0.20-38    2018-11-04 [2] CRAN (R 3.6.0)
#>  lazyeval        0.2.2      2019-03-15 [1] CRAN (R 3.6.0)
#>  listenv         0.7.0      2018-01-21 [1] CRAN (R 3.6.0)
#>  lubridate     * 1.7.4      2018-04-11 [1] CRAN (R 3.6.0)
#>  magrittr        1.5        2014-11-22 [1] CRAN (R 3.6.0)
#>  Matrix          1.2-17     2019-03-22 [2] CRAN (R 3.6.0)
#>  memoise         1.1.0      2017-04-21 [1] CRAN (R 3.6.0)
#>  mgcv            1.8-28     2019-03-21 [2] CRAN (R 3.6.0)
#>  munsell         0.5.0      2018-06-12 [1] CRAN (R 3.6.0)
#>  nlme            3.1-139    2019-04-09 [2] CRAN (R 3.6.0)
#>  pillar          1.4.1      2019-05-28 [1] CRAN (R 3.6.0)
#>  pkgbuild        1.0.3      2019-03-20 [1] CRAN (R 3.6.0)
#>  pkgconfig       2.0.2      2018-08-16 [1] CRAN (R 3.6.0)
#>  pkgload         1.0.2      2018-10-29 [1] CRAN (R 3.6.0)
#>  plyr            1.8.4      2016-06-08 [1] CRAN (R 3.6.0)
#>  prettyunits     1.0.2      2015-07-13 [1] CRAN (R 3.6.0)
#>  processx        3.3.1      2019-05-08 [1] CRAN (R 3.6.0)
#>  ps              1.3.0      2018-12-21 [1] CRAN (R 3.6.0)
#>  purrr         * 0.3.2      2019-03-15 [1] CRAN (R 3.6.0)
#>  R6              2.4.0      2019-02-14 [1] CRAN (R 3.6.0)
#>  Rcpp            1.0.1      2019-03-17 [1] CRAN (R 3.6.0)
#>  readr         * 1.3.1      2018-12-21 [1] CRAN (R 3.6.0)
#>  remotes         2.1.0      2019-06-24 [1] CRAN (R 3.6.0)
#>  reshape2        1.4.3      2017-12-11 [1] CRAN (R 3.6.0)
#>  rlang           0.4.0      2019-06-25 [1] CRAN (R 3.6.0)
#>  rmarkdown       1.13       2019-05-22 [1] CRAN (R 3.6.0)
#>  rprojroot       1.3-2      2018-01-03 [1] CRAN (R 3.6.0)
#>  scales          1.0.0      2018-08-09 [1] CRAN (R 3.6.0)
#>  sessioninfo     1.1.1      2018-11-05 [1] CRAN (R 3.6.0)
#>  storr           1.2.1      2018-10-18 [1] CRAN (R 3.6.0)
#>  stringi         1.4.3      2019-03-12 [1] CRAN (R 3.6.0)
#>  stringr         1.4.0      2019-02-10 [1] CRAN (R 3.6.0)
#>  testthat        2.1.1      2019-04-23 [1] CRAN (R 3.6.0)
#>  tibble          2.1.3      2019-06-06 [1] CRAN (R 3.6.0)
#>  tidyr         * 0.8.3      2019-03-01 [1] CRAN (R 3.6.0)
#>  tidyselect      0.2.5      2018-10-11 [1] CRAN (R 3.6.0)
#>  usethis         1.5.0      2019-04-07 [1] CRAN (R 3.6.0)
#>  withr           2.1.2      2018-03-15 [1] CRAN (R 3.6.0)
#>  xfun            0.8        2019-06-25 [1] CRAN (R 3.6.0)
#>  yaml            2.2.0      2018-07-25 [1] CRAN (R 3.6.0)
#> 
#> [1] /Users/shik544/R
#> [2] /usr/local/Cellar/r/3.6.0_3/lib/R/library

#> Local:    master /Users/shik544/Box Sync/Projects/forte_project/fortebaseline
#> Remote:   master @ origin (git@github.com:ashiklom/fortebaseline.git)
#> Head:     [7dc6465] 2019-07-01: Add results text about lai x pft figure